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Abstract 

We report studies of ultra-high energy cosmic ray composition via analysis of depth of airshower 
maximum (X max ), for airshower events collected by the High Resolution Fly's Eye (HiRes) ob- 
servatory. The HiRes data are consistent with a constant elongation rate d <X m ax> I 'd(log(E)) 
of 47.9 ± 6.0(stat.) ± 3.2(syst.) g/cm 2 /decade for energies between 1.6 EeV and 63 EeV, and are 
consistent with a predominantly protonic composition of cosmic rays when interpreted via the 
QGSJET01 and QGSJET-II high-energy hadronic interaction models. These measurements con- 
strain models in which the galactic-to-extragalactic transition is the cause of the energy spectrum 
"ankle" at 4 x 10 18 eV. 

PACS numbers: 98.70.Sa, 95.85.Ry, 96.50.sb, 96.50.sd 
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The observation of a break in the cosmic ray energy spectrum at approximately 6 x 
10 19 eV jl-3| provides evidence that the highest energy cosmic rays are both extragalactic 
and protonic 4j, |5|. Direct evidence for a proton-dominated composition from airshower 
data would lend further support to this model, as would the observation of a transition from 
(heavy) galactic to (light) extragalactic cosmic rays at lower energies. A second feature, the 
"ankle" of the energy spectrum at 4 x 10 18 eV may be indicative of this transition or it 
may further strengthen the model in which the end of the cosmic ray spectrum is shaped 
by interactions with the cosmic microwave background [6] . Composition studies can provide 
decisive evidence in the choice between interpretations. 

An important clue to chemical composition which is accessible to air fluorescence obser- 
vatories is the depth of shower maximum X max of cosmic ray induced extensive airshowers. 
Simple arguments 3, 0] show that the average value of airshower maximum <X max > will 
depend logarithmically on the primary energy and atomic mass, and that the elongation rate 
d <X max > /dlogE will be constant for unchanging primary compositions. Further, to first 
order, a nucleus-induced shower may be thought of as a superposition of showers induced 
by single nucleons. Therefore due to averaging effects we also expect the width of the X max 
distribution at a given energy to be sensitive to the atomic mass of the primary. 

The two fluorescence observatories of the High-Resolution Fly's Eye collected data in 
stereoscopic mode between December 1999 and April 2006. Located on the U.S. Army 
Dugway Proving Ground in Utah, at a mean elevation of 1,575 m MSL, a mean latitude of 
40.16° N, and separated by 12.6 km, the observatories operated on clear moonless nights. 
Each detector consisted of an array of telescope modules, each module included a mirror of 
3.7 m 2 effective area which focused UV light from airshowers on a 16 x 16 photomultiplier 
tube (PMT) camera. The field of view of each PMT subtended a one degree diameter cone 
of the sky. The HiRes-I detector covered nearly 360° in azimuth, 3°-17° in elevation and 
was read out by means of sample-and-hold electronics, while the HiRes-II detector covered 
3°-31° in elevation and was read out by a custom FADC system J9]. 

The calibration of the HiRes telescope modules has been described previously [l0|. A 
portable Xenon flash lamp (~ 0.5% stability) was used to illuminate each mirror monthly. 
Between Xenon runs, nightly calibrations were performed using YAG laser light delivered 
to the cameras via optical fiber A pulsed nitrogen laser was fired into the atmosphere 
from various locations within 3.5 km of the two detector sites. An overall accuracy of ~ 10% 
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RMS is achieved in the HiRes photometric scale. 

Steerable lasers fired patterns of shots which covered the aperture of the HiRes fluores- 
cence detectors, in order to monitor UV attenuation in the atmosphere. The vertical aerosol 
optical depth (VAOD) was measured to be 0.04(mean) ±0.02(RMS) fl^l. corresponding 
to a mean correction of ~15% upward in energy for an event 25 km distant from the obser- 
vatory. In the present analysis, the steerable laser measurements were used to compile an 
hourly database of the atmospheric parameters. 

A mirror trigger was initiated if a sufficient number of PMTs were in temporal and spatial 
coincidence, then a stereo data set was obtained by the time-matching of HiRes-I and HiRes- 
II triggers. Geometrical reconstruction of stereo events proceeded by determination of the 
shower-detector plane from each HiRes site, the intersection of these two planes was taken 
to be the shower core trajectory. The resolution in the shower zenith angle is 0.6°, and the 
resolution in R p (distance of closest approach to the detector) is 1.2%. 

Hit information from multiple tubes in the HiRes-II data only are sorted into discrete 
time bins. In each bin, FADC signals are converted to a number of photoelectrons N pe , 
then adjusted for the effective area of each bin as determined by ray tracing. The geometry 
of the shower and the atmospheric databases are used to determine the atmospheric slant 
depth X for each shower bin. Shower segments that have emission angles within 5° of a 
bin's pointing direction are flagged and not used for fitting due to excessive Cherenkov light 
contamination. The N pc profile is then converted to a profile of the fluorescence light at 
the shower by a routine that simulates the light production and propagation through the 
atmosphere, and subtracts the Cherenkov contribution. 

The intensity of fluorescence light emitted from an airshower is pro portional to the total 
ionization energy deposited by the charged particles in the shower 141 ] . We use the average 
fluorescence yield of several groups [l5-17] and the spectral distribution given by Ref. 18], 
along with the average dE/dX determined from CORSIKA simulations 19|] to determine 
the number of charged particles in the shower as a function of slant depth. 

The shower profile is then fit to a Gaussian function of the age parameter s{X) = 3X/ (X+ 
2A max ) in order to determine the airshower energy and X max . (Alternatively fitting by the 
Gaisser-Hillas parametrization |20[ had little overall effect on the analyses and conclusions 
presented here.) Further details of the reconstruction used in this analysis are contained in 



Ref. 



2l|. 
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Ultra-high energy cosmic ray composition studies require a detailed comparison of data 
with the predictions of airshower models and are hence model dependent. In practice, these 
models are airshower Monte Carlo simulations. In order to completely understand the effects 
of the geometrical aperture of the detector, HiRes applies a full detector simulation to the 
Monte Carlo events. 



HiRes uses libraries of simulated proton- and iron-induced airshowers generated 



3y the 



CORSIKA_6.003 (6.501) [22j package, using the QGSJET01 [23j, QGSJET-II [24|, and 



SIBYLL 
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26j hadronic interaction models and the EGS4 [27| electromagnetic interaction 



driver. The number of particles versus slant depth is recorded at 205 points along each 
shower. 

Detector simulation proceeds by drawing an event from the shower library, assigning it 
a random core location, zenith and azimuthal angle and determining if it can trigger the 
detector. The number of charged particles at many discrete points along the shower is de- 
termined, and fluorescence light is then propagated from the shower to the detector. Light 
attenuation by the atmosphere is realistically simulated by using an hourly database de- 
scribing the measured aerosol content, temperature, and pressure. Ray tracing is performed 
to determine which phototubes see the light, allowing for photocathode response and in- 
active space between PMTs. An electronics simulation then determines the pulse height 
and time for each tube for HiRes-I and forms a FADC waveform for HiRes-II. For all tubes 
the channel gains, DAC settings, thresholds and channel variances are simulated by using 
hourly database information. The same trigger algorithms used in hardware are simulated, 
and if either detector would have been triggered the simulated data is written to disk in the 
identical format as real data, allowing the study of Monte Carlo events by the same analysis 
chain. The Monte Carlo set for each hadronic interaction model contains approximately 
20 times the number of reconstructed stereo events as the data. 

The major challenge in studying cosmic ray composition by the <X max > technique lies in 
understanding the systematic biases that occur during reconstruction and event selection. 
Low-energy showers which are nearby the detector may reach X max above the field of view of 
the mirrors, thus biasing a data set towards deeper showers. High-energy showers with small 
zenith angles are likely to reach maximum below the field of view of the mirrors, resulting 
in a bias towards shallow showers. 

It is useful to divide the types of biases which can occur into two types, called "accep- 
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FIG. 1: Difference between HiRes-II (Xjj) and HiRes-I (Xj) X max for HiRes stereo data (points) 
overlaid with QGSJET-II proton Monte Carlo. X max bracketing by HiRes-I is required. 

tance biases" and "reconstruction biases". Acceptance biases are due to events which fail 
reconstruction altogether, including detector triggering and event selection effects. Recon- 
struction biases are due to events which are successfully reconstructed, but with the wrong 
Xmax- The strategy in the following analysis is to choose event selection cuts which minimize 
the reconstruction bias, and make the acceptance bias as independent of energy as possible. 

After geometrical reconstruction and fitting of the shower profiles to obtain energy and 
Xmax, a set of cuts are applied in order to select an appropriate event sample. The chance 
probability that the event is due to noise must be less than 1%, and the x 2 /DOF of the 
fit must be less than 4. Data must have been taken in good weather conditions. Fit 
uncertainty in the zenith angle must be less than 2°, the fit uncertainty in X max must be 
less than 40 g/cm 2 , and the angular RMS with respect to the event plane of hit PMTs must 
be greater than 0.15°. The zenith angle of the event itself must be less than 70°, and the 
R p with respect to HiRes-II must be at least 10 km. Events are required to have X max 
bracketed by the HiRes-II field of view, and have a shower-detector plane angle between 
40° and 130°. Finally, events with energies 18.2 < logf^-E^eV)) < 19.8 are selected for this 
analysis. A total of 815 events pass all cuts. 

In Fig. [1] the resolution in X max of data and Monte Carlo events are compared by plotting 
the difference between X max as measured by HiRes-I and HiRes-II. The agreement is excel- 
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lent, including the asymmetry caused by HiRes-I covering only half the range in elevation 
angle. This supports the use of Monte Carlo to determine X max resolution, which is found 
to be better than 25 g/cm 2 over most of the HiRes energy range. 
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FIG. 2: Top: X max overlay of HiRes data (points) with QGSJET-II proton Monte Carlo airshowers 
after full detector simulation. Bottom: X max overlay of HiRes data (points) with QGSJET-II iron 
Monte Carlo airshowers. 

After application of the cuts above, we compare the distribution in X max to the pre- 
dictions obtained from simulated proton and iron showers. Fig. [2] shows the excellent 
overall agreement between the HiRes stereo data and the predictions of the QGSJET-II 
proton Monte Carlo. In Fig. [3l we compare the development with energy of the mean 
X max in HiRes data with the predictions (after full detector simulation) for QGSJET01, 
QGSJET-II and SIBYLL proton and iron. The data agrees best with the QGSJET-II pro- 
ton prediction, with a x 2 = 6.9/8 df. In a linear fit with \ 2 = 5.2/(6 df) we measure 
47.9 ± 6.0(stat.) g/cm 2 /decade for the elongation rate. Systematic effects are considered 
below. 
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FIG. 3: HiRcs stereo <X max > compared with the predictions for QGSJET01, QGSJET-II and 
SIBYLL protons and iron after full detector simulation. The number of events in each energy bin 
is displayed below the data point. 

Due to detector, reconstruction, and event selection acceptance effects, the proton and 
iron "rails" in Fig [3] are shifted relative to the raw CORSIKA predictions. However, we find 
that the shift in mean X max for QGS JET01 and QGSJET-II protons is independent of energy 
to approximately 1.8 g/cm 2 /decade. We assign a systematic uncertainty of 2.7 g/cm 2 /decade 
to the elongation rate based on small variations of event selection cuts. Uncertainties in 
the energy do not have a large effect on elongation rate results due to the logarithmic 
energy scale. The choice of VAOD was the main systematic in a previous elongation rate 



analysis 28j], however the use of an hourly atmospheric database in the present analysis 



renders this source of systematics negligible. 



The phototube pointing directions have been confirmed by studies using stars 29j to 
within 0.3°, corresponding to a shift in X max of approximately 15 g/cm 2 . Averaging over 
mirrors, this contributes a net uncertainty of 3.3 g/cm 2 to the value of <X max >. The 
subtraction of the Cherenkov light from the phototube signal can introduce an uncertainty 



in X max due to uncertainties in electron multiple scattering. Previous studies [28| in which 



the width of the Cherenkov beam was varied by 2° (1 a) indicated negligible effect on the 



8 



elongation rate or absolute value of <X max >. Finally, a systematic uncertainty of 0.7 g/cm 2 
is assigned to the absolute value of <X max > in the predictions due to Monte Carlo statistics. 

The fluctuations of X max as a function of energy are also a probe of primary particle com- 
position. Because the distributions tend to be both asymmetric and possess non-Gaussian 
tails, care must be taken to use a suitable definition of the X max width. The uncorrected 
RMS and sample standard deviations are biased estimators of the width [30] and tend to be 
subject to large fluctuations in distributions with broad tails. 

In order to focus attention on the center of the X max distribution and reduce sensitivity 
to fluctuations in the tails, the width is quantified as the width ax of a unbinned likelihood 
fit to a Gaussian of a distribution truncated at 2 x RMS. The results of this analysis 
applied to both the HiRes data and to QGSJET-II proton and iron Monte Carlo are shown 
in Fig. HI The HiRes X max width data are consistent with the predictions of QGSJET-II 
protons. The width of the X max distribution of protons within the QGSJET01 model tends 
to be about 5 g/cm 2 broader than that of QGSJET-II, while SIBYLL protons are 2-3 g/cm 2 
narrower than those of QGSJET-II. Both of these shifts are small compared with statistical 
uncertainties, particularly at the highest energies. 

In summary, the HiRes data are consistent with a constant elongation rate of 47.9 ± 
6.0(stat.)±3.2(syst.) g/cm 2 /decade above 1.6 EeV, and thus with an unchanging composition 
across the ankle. This places strong constraints on models in which the ankle is the result 
of a transition from heavy galactic to light extragalactic cosmic rays. 

Of the hadronic interaction models tested, the best agreement is with the QGSJET-II 
pure proton model. Within current uncertainties the data are completely consistent with 
this model, and close to QGSJET01 pure protons. Comparison with SIBYLL suggests a 
mixture dominated by light elements. The observed constant elongation rate would imply 
that this mixture is unchanging, or at most steadily changing over nearly two orders of 
magnitude spanning the energy spectrum ankle. 

The present analysis, taken together with the HiRes spectral results jl, 3] on the shape 
and location of the high-energy cutoff and ankle, suggests the simple picture in which cosmic 
rays above 1 EeV are protons of extragalactic origin and the end of the energy spectrum is 
shaped by interactions with the cosmic microwave background. 

PHY-9904048, PHY-9974537, PHY-0098826, PHY-0140688, PHY-0245428, PHY- 
0305516, PHY-0307098, and by the DOE grant FG03-92ER40732.We gratefully acknowl- 
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FIG. 4: Results of fitting HiRes stereo data X max distribution to Gaussian truncated at 2 x RMS 
(black points). Superimposed are expectations based on QGSJET-II proton (squares) and iron 
(triangles) Monte Carlo. Monte Carlo points are shown with small offsets in energy to provide 
separation. 
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Note added. - In a recent paper [3l|], the Pierre Auger collaboration draws different 
conclusions as to the composition of the highest energy cosmic rays. 
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